Dynamical heterogeneity in a model for permanent gels: Different behavior of 

dynamical susceptibilities 
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We present a systematic study of dynamical heterogeneity in a model for permanent gels, upon 
approaching the gelation threshold. We find that the fluctuations of the self intermediate scattering 
function are increasing functions of time, reaching a plateau whose value, at large length scales, 
coincides with the mean cluster size and diverges at the percolation threshold. Another measure 
of dynamical heterogeneities, i.e. the fluctuations of the self-overlap, displays instead a peak and 
decays to zero at long times. The peak, however, also scales as the mean cluster size. Arguments 
are given for this difference in the long time behavior. We also find that non-Gaussian parameter 
reaches a plateau in the long time limit. The value of the plateau of the non-Gaussian parameter, 
which is connected to the fluctuations of diffusivity of clusters, increases with the volume fraction 
and remains flnite at percolation threshold. 

PACS numbers: 82.70.Dd, 64.60.Ak, 82.70.Gg 



I. INTRODUCTION 

In the context of the glass transition the concept of dy- 
namical heterogeneities has been very fecund [1-12]. In 
glassy systems the correlated motion of particles mani- 
fest as significant fluctuations around the average dynam- 
ics, strongly increasing as the transition is approached. 
These heterogeneities in the dynamics have been studied 
quantitatively via the so-called dynamical susceptibility 
[3], Xi{t) ^ N{(F^{t))-{F{t))^), obtained as the fluctua- 
tions of a suitable time dependent correlator F{t) (where 
N is the number of particles and (...) is the ensemble 
average). Two quantities are usually considered: The 
fluctuations of the self intermediate scattering functions 
(ISF) [6, 13] X4(fc,t) = N[{\^,ik,t)\^}-{^,ik,t)y] 
usually measured in numerical simulations, and the fluc- 
tuations of the time dependent overlap [14-16] xi{t) = 
N[{q{t)'^) — {q{t))'^], which, first introduced inp-spin glass 
models [3], has been calculated also within mode cou- 
pling theory [3, 13]. For the fluctuations of the overlap, 
the role of the inverse of the wave vector k is played 
by the parameter a characterizing the overlap function, 
which is different from zero only if a particle has moved 
a distance less than the fixed value a. In usual glassy 
systems the behavior observed in the dynamical suscep- 
tibility is essentially the same despite of different choices 
of F{t) [4, 6, 13]: X4(0 grows as a function of the time, 
reaches a maximum and then decreases to a constant, 
consistently with the transient nature of the dynamical 
heterogeneities. Some differences in the k dependence of 
these two quantities were however found in a model for 
glasses [17]. 

Recently, dynamical heterogeneities have been stud- 
ied in other complex systems, such as granular media 



[18-21] and attractive colloidal systems [22-26], where 
behaviors qualitatively similar to that found in glasses 
are observed. In particular, Ref.[26] reports a systematic 
study of the dynamic susceptibility in colloidal systems 
along the attractive glassy line. Typically the dynami- 
cal susceptibility, defined as the fluctuations of the self 
ISF, displays a well pronounced peak. However, in the 
attraction-dominated limit, the dependence on both time 
and wave vector markedly differs from that in standard 
repulsion-dominated systems (hard-sphere limit). 

In a recent letter [27] we have studied the behavior of 
the dynamical susceptibility, X4, dcflned as the fluctua- 
tions of the self ISF, in a model for permanent gel, where 
bonds are modeled using a finitely extendable non linear 
elastic (FENE) potential [29, 30] between neighboring 
particles. It was found that the behaviour of X4(A:,<) is 
drastically different from that found in glasses. In fact 
it grows in time until it reaches a plateau in the limit 
of large time i, without decaying to 1. The value of the 
plateau in the limit of low wave vector, /c — > 0, was in 
fact found to coincide with the mean cluster size. As a 
consequence, as the system approaches the gel transition 
(i.e. the percolation threshold), the value of the plateau 
diverges. For a fixed value of fc, the value of the plateau 
coincides with the mean cluster size up to a linear size of 
the order of the inverse of k. Therefore, for any fc > 
{k > kmin in our study), the plateau never diverges: it 
decreases as k increases and eventually goes to one. 

Here we present a systematic study of this FENE 
model for permanent gels [27]. Moreover, we compare 
the behavior of the fluctuations of the self ISF and of 
the self-overlap, and flnd a marked difference between 
the two ones. The first one, as mentioned above, is an 
increasing function of the time and tends to a plateau. 
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whereas the second one reaches a maximum and then de- 
creases. However, the value of the maximum scales as the 
value of the plateau of the fluctuations of the self ISF, 
with the same critical exponent of the mean cluster size. 

The reason why these two quantities differ so drasti- 
cally in the long time limit is the following: the fluctua- 
tion of the overlap is related to the correlations between 
the event that a monomer has moved a distance less then 
a in a time interval t and the event that another monomer 
has also moved a distance less then a in the same interval 
t. In the long time t all particles have moved a distance 
larger than a therefore such correlations are zero. On 
the other hand in the long time limit the fluctuation of 
the self ISF is related to the correlation of the distance 
separating monomers i and j at time and the distance 
between the same monomers at time t. This quantity is 
different from zero if the particles i and j are in the same 
cluster. 

Although the long time limit of the two quantities 
XA{k,t) and x?(a;i) is different in gels with permanent 
bonds, they have in common not only the property that 
plateau and maximum scale in the same way, but also 
one key feature which is the strong length scale depen- 
dence: The peak of the fluctuations of the self part of 
the overlap and the plateau of the fluctuations of the self 
ISF decreases strongly as the wave vector k (or 1 /a) in- 
creases, which is the sign that clusters of bonded particle 
dominate the dynamics. The same feature is also valid 
for strong colloidal gels [31]. This strong length scale de- 
pendence of the dynamical susceptibility seems to be the 
distinctive sign of permanent or strong colloidal gelation, 
compared with the (attractive or repulsive) glass transi- 
tion. 

Finally, we measure the non-Gaussian parameter a2 
and find that, due to the presence of clusters, it is dif- 
ferent from zero also in the long time limit. However, its 
plateau value, which is connected to diffusivity, remains 
finite upon approaching the transition. 

In Sect. II we introduce the model used and give the 

details of the numerical simulations. We analyze the self 
ISF and its fluctuations in Sect. HI and the self-overlap 
and its fluctuations in Sect. IV. The mean squared dis- 
placement and the non- Gaussian parameter are discussed 
in Sect.V, whereas the Sect. VI contains the conclud- 
ing remarks. Finally, in Appendix A we investigate the 
static properties of the sol-gel transition, corresponding 
to the percolation of permanent bonds between particles 
[32, 33]. 



II. MODEL AND NUMERICAL SIMULATIONS 

We consider a 3d system of A'' particles interacting 
with a soft potential given by Weeks-Chandler-Andersen 



(WCA) potential [28]: 



rrWCA _ r 4e[(a/ry)'' - (^/ry)« + \], nj < 2^^ 
""'^ - \ 0, nj > 21/6^ 

(1) 

where is the distance between the particles i and j. 

After the equilibration, particles distant less than Rq 
are linked by adding an attractive potential: 



jjFENE ^ f -0.5fco-Roln[l - inj/Ro)^], nj < Rq 

\ oo, > i?o 

(2) 

representing a finitely extendable nonlinear elastic 
(FENE). The FENE potential was firstly introduced in 
Ref. [29] and is widely used to study linear polymers [30] . 
We choose /cq = 30e/a-^ and i?o = l-5cr as in Ref. [30] in 
order to avoid any bond crossing and to use an integra- 
tion time step At not too small [34]. The introduction of 
the FENE potential leads to the formation of permanent 
bonds among all the particles whose distance at that time 
is smaller than Rq. 

We have performed moleciilar dynamics simulations of 
this model: The equations of motion were solved in the 
canonical ensemble (with a Nose-Hoover thermostat) us- 
ing the velocity- Verlct algorithm [35] with a time step 
At = 0.001(5r, where 5t ~ aimje)^^'^ is the standard 
unit time for a Lennard- Jones fluid and m is the mass 
of particle. We use reduced units where the unit length 
is cr, the unit energy is e and the Boltzmann constant 
]tB is set equal to 1. We use periodic boundary condi- 
tions, and average all the investigated quantities over 32 
independent configurations of the system. 

The temperature is fixed at T = 2 and the volume 

fraction <^ ~ ira^N/QL^ (where L is the linear size of the 
simulation box in units of a) is varied from (j) = 0.02 
to <f) = 0.12. Using the percolation approach, we iden- 
tify the gel phase as the state in which there is a per- 
colating cluster [32, 33]. A finite size scaling analysis 
is presented in the Appendix, showing that this transi- 
tion is in the universality class of random percolation. 
We find that the threshold is = 0.09 ± 0.01. In 
particular, we obtain that the cluster size distribution, 
Us ~ s~'^ for <j> — (j>c with r = 2.1 ± 0.2, the mean 
cluster size S{(j)) = J2^'^''^s/J2^''^s ^ {4>c — <t>)~'^ with 
7 = 1.8±0.1, and the connectedness length ^ ^ ((Ac^ </>) 
with v = 0.88 ± 0.01. In the following we fix the number 
of particles, N = 1000. 

Due to the introduction of bonds, spatial correlations 
appear at low wave vectors. Although these correlations 
increase as a function of the volume fraction, the low k 
limit of the static structure factor S{k) is always small 
compared to the number of particles, and no phase sep- 
aration is observed. 
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III. SELF INTERMEDIATE SCATTERING 
FUNCTION AND ITS FLUCTUATIONS 

Relevant information on the relaxation dynamics over 
different length scales can be obtained from the self In- 
termediate Scattering Functions (ISF) Fs{k,t): 



Fs{k,t) = [($,(fc,t))] 



(3) 



where (...) is the thermal average over a fixed bond con- 
figuration, [...] is the average over 32 independent bond 
configurations of the system, and 



N 



,ik-{ri{t)-ri{0)) 



(4) 



In Fig.l, Fs{k, t) is plotted as a function of t for differ- 
ent 4>, respectively for kmin = 27r/L ~ 0.35 (main frame) 
and fc ~ 7 (inset). At the smallest wave vector, for very 
low values of the volume fraction, the self ISF decays to 
zero following an exponential behavior. As the volume 
fraction is increased towards the percolation threshold, 
we observe the onset of a stretched exponential decay, 
g-(*/T) ^ with P decreasing as a function of cj) (for in- 
stance P = 0.75 ± 0.01 for = 0.07 and P = 0.58 ± 0.02 
for (j) = 0.085). The cluster size distribution has started 
to widen towards the percolation regime (see Appendix) , 
and therefore, over sufficiently large length scales, the 
behavior of Fs{k,t) is due to the contribution of differ- 
ent relaxation processes, characterized by different relax- 
ation times, whose superposition produces a detectable 
deviation from an exponential law. Near the transition 
threshold the long time decay is characterized by a power 
law behavior, indicating that the relaxation over this 
length scale is controlled by the formation of the perco- 
lating cluster, with a critically growing relaxation time. 
If the volume fraction increases further, the decay be- 
comes slower and slower, showing a logarithmic behav- 
ior. These features of the dynamics well reproduce the 
experimental findings [36] . Moreover they agree with re- 
sults obtained via numerical simulations of different gela- 
tion models [37-39]. At large wave vectors (see Inset of 
Fig.l) and low volume fractions, Fs{k,t) decays to zero 
as e"^*/"^^ (continuous curves in figure), corresponding 
to the ballistic regime of particle motion. 

From the self ISF we calculate the structural relaxation 
time, T{k, (j)), defined as the time for which Fs{k, T{k)) ~ 
0.1. In Fig. 2, T{k,(j)) is plotted for different values of k 
as a function of the volume fraction (f). For k = kmim we 
find that T{kmim 4') is well fitted by a power law diverg- 
ing at the gelation threshold with an exponent / ^ 1.22. 
Increasing k, no divergence is observed at the thresh- 
old, signalling that no structural arrest occurs over length 
scales less than the box size, L. 

In Fig. 3 and in its inset we plot respectively k'^r^k, (j)) 
and kT{k, (f>) as a function of the wave vector for different 
volume fractions. The inset of Fig. 3 shows that r ~ l//c 
for large wave vectors, refiecting the ballistic diffusion 




FIG. 1: (Color online) Main frame: Self ISF for cl> = 0.02, 
0.07, 0.08, 0.09, 0.1 (from bottom to top) and k ~ 0.35 as a 
function of time t. The lines are fitting curves: For (f) « (j)c 
the decay is well fitted by an exponential behavior (dashed 
line); if cf> approaches to (j)c a stretched exponential decay 
appears, with /3 = 0.75 ± 0.01 for = 0.07 (dotted hne). For 
4> ~ 0.09 the decay is well fitted by a power law ~ t~'^ with 
c = 0.65 ± 0.03 (full line). Inset: Self ISF for A: ~ 7 and the 
same volume fractions of main frame. 



for short times (see Sect.V). Interestingly in the limit of 
small wave vectors fc^r does not tend to a constant. This 
unusual result is essentially due to the fact that the the 
non-Gaussian parameter [40], a2(t) = (^^^yi — 1, is 
non zero in the long time limit, as discussed in details in 
Sect.V. In this case the Gaussian approximation of the 
probability distribution of particle displacements is not 
valid, and the self ISF Fs{k,t) cannot be written as a 
Gaussian even in the limit of small wave vector. 

We now analyze and discuss the behaviour of the fiuc- 
tuations of the self ISF, i.e. the dynamical susceptibility: 



X4{k,t)^N (|$,(A;,i)|2)-(ci>,(fc,t))2 



(5) 



In Fig. 4 Xii^Tt) is plotted for k = kmin and different 
volume fractions. Differently from the behavior typically 
observed in glassy systems, we find that, for cf) < (j^c, 
Xiik, t) is a monotonically increasing function of the time 
tending to a plateau in a time of the order of the relax- 
ation time T(kmin)- The value of the plateau diverges 
as the mean cluster size as the percolation threshold is 
approached [27]. For (j) > (j)c the system is out of equilib- 
rium, X4(/c,i) continues increasing as a function of time, 
without reaching any asymptotic value within the simu- 
lation time. We briefiy discuss the main arguments ex- 
plaining the above result, presented in Ref [27], where 
it was in fact shown that, for k —^ and t ^ oo, the 
dynamical susceptibility X4(fc,t) tends to the mean clus- 



ter size. We define Xas{k,(f)) = lim 



N- 



, lim* 



.X4{k,t). 
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FIG. 2: (Color online) Structural relaxation time r(fc, (j)) as 
a function of the volume fraction, for wave vector k ~ 0.35, 
0.6, 1.0, 2.0, 3.0, 7.0 (from top to bottom). The fuU line is 
the fitting curve: T{kmin,4>) ~ {4>c — 0)"^, with / ~ 1.22. 
Dashed lines are eye guides. 



FIG. 4: (Color online) Dynamical susceptibility, Xiik,t), as a 
function of time for k = kmin and different volume fractions 
(?!) = 0.05, 0.06, 0.07, 0.08, 0.09 (from bottom to top). 
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FIG. 3: (Color online) Main frame: k'^T{k,(j)) as a function 
of k, for (l> = 0.05, 0.06, 0.07, 0.08, 0.09, 0.1 (from bottom to 
top). Inset: fcr(A;, 0) as a function of k for the same volume 
fractions of main frame. 



Being limt^oo {^s{k,t)) = 0, we have 



Xas{k,(j)) = lim — 



N 

E 



Cij (k) 



(6) 



where Cy(fc) ^ linit^^(e*^-(*''(*)-*'^(*))e-*^-(^'(°)-*'^"(°») = 
j^gi/c (ri-rj)^|2^ Here we have used the fact that, for large 
enough time t, the term e"''''^''"^*-'"''^*^'^-' is statistically 



independent from e^''^'(''*(*'^^''j*^°)\ so that we can 
factorize the thermal average. We separate the sum over 
connected pairs (7^ = 1, i.e. pairs belonging to the 
same cluster) and disconnected pairs (7^ — 0, i.e. pairs 
belonging to different clusters), so that: 



Xas(fc,0) 



lim — 



(7) 



N 



lijCij{k) 



1 

TV 



N 



E (l-7.,)a,(fc) 



If particles i and j are not connected, for any fixed value 
of A; > 0, the quantity C„ (fc) is O^i/N"^) [41]. As there 
are at most iV^ disconnected pairs, the second term of 
the r.h.s. of Eq.(7) is 0{1/N), and can be neglected in 
the thermodynamical limit. 

For (j) < (pc, clusters will have at most a linear size of 
order ^, so that the relative distance \fi 
nected particles will be smaller than ^. 
have limfc_»o7y Cy (fc) = 7^^- and 



fj I of con- 
Therefore we 



lim Xas{k,4>) 



1™ TT 



N 



E 



(8) 



where S is the mean cluster size. As shown in Ref.[27], 
numerical data confirm this result. 

In Fig. 5 Xi{k,(i>) is plotted for = 0.09 and differ- 
ent wave vectors. For each value of the wave vector, 
X4{k,(j)) reaches a plateau after a characteristic time of 
the order of the relaxation time T(fc). The asymptotic 
value Xasik, (f>c) at low wave vectors follows a scaling be- 
havior as a function of k (Inset of Fig. 5): at the tran- 
sition threshold the exponent is 2.03 ± 0.02, consistent, 
within the numerical accuracy, with the prediction 2 — t] 
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FIG. 5: (Color online) Main Frame: Dynamical susceptibil- 
ity, Xiikjt), as a function of time for (j) = 0.09 and k = 0.35, 
0.61, 0.99, 1.40, 2.10, 3.96 (from top to bottom). Inset: 
Asymptotic values of the susceptibility, Xas{k,(l)c) as a func- 
tion of the wave vector k. Data are fitted with a power law 
~ /j-2 03±o.02^ agreement with the exponent 2 — r; of random 
percolation. 



of random percolation [42]. This result shows that if 
one varies the wave vector k (and 27r/fc > a) the dy- 
namical susceptibility is able to detect the self-similarity 
of the structure of the system due to the percolation 
transition. Using scaling arguments [27], we can write 
Xas{k,(t)) = k^~^f{k£^) where f{z) is a function, which 
tends to a constant for small z, whereas it behaves as z''^'^ 
for large values of z. As shown in Ref.[27], data support 
this scenario. All these results coherently show how in 
the present system the asymptotic value of the dynamical 
susceptibility can be related to the cluster size. Not only 
our results indicate that the percolation exponents can 
be measured in a direct way, by developing techniques to 
measure the dynamical susceptibility, but they also state 
that the asymptotic value of the dynamical susceptibility 
plays the same role as the static scattering function near 
the liquid gas critical point. 



IV. SELF OVERLAP AND ITS FLUCTUATIONS 

In the context of glassy systems, a time-dependent or- 
der parameter was introduced [14-16], which measures 
the number of "overlapping" particles in two configura- 
tions separated by a time interval t, 

q{t) = ]^ y dfidf2p{ri,Q)p{r2,t)w{\fi - r2\) 

i 3 



where p(r, t) = S{f— ri{t)) is the density in r at time 
t, and 'w{\fi — r2|) is an "overlap" function that is 1 for 
— ^2 1 < 1 and zero otherwise [43] . 
In Ref.[4] the authors separate q into self and distinct 
components, q{t) = qs{t) + qoit)- The self part is given 
by: 

isit)-^J2''i\m^nm (10) 

i 

which corresponds to terms of Eq.(9) with i = j, and 
measures the number of particles that move less than a 
distance a in a time interval t. In Ref.[4] it was shown 
that on average the dominant term is given by the self 
part. 

Here we measure Qs{t) — [{Qs{t))] for two choices of a, 
0.15 and 3, respectively corresponding to 1/a >> kmin 
and 1/a ^ kmin- Qs{t) is plotted in Fig. 6 for different 
values of the volume fraction. We have verified by nu- 
merical calculations that for small enough a, the relevant 
contribution to Q{t) — [{q{t))] is given by Qs{t), since 
the probability that a particle replaces within a radius a 
another particle is small. For all the values of a and of 
considered, Qs{t) at long times is well fitted by a power 
law. 

Another interesting method to investigate the spatially 
heterogeneous dynamics, generally used in glassy sys- 
tems, is the measure of the dynamical susceptibility ob- 
tained by the fluctuations of the time dependent over- 
lap [14-16] x?(a,i) = N[{qit)^) - {q{t)y], where q{t) is 
given by Eq.(9). In glassy systems this quantity essen- 
tially presents the same features as the fluctuations of 
the self ISF. 

Here we measure the fluctuations of the self part of the 
overlap: 

x'^ia,t)=N[{qs{tf)-{qsit)n (H) 

for different choices of a, ranging from 0.15 to 3. In Fig. 7 
we plot Xi {(^7 for a = 3 and different values of 0. 

We see that differently from the fluctuations of the self 
ISF, here x^(a, t) displays a peak, whose value increases 
and diverges as the the gel transition is approached. In- 
deed, the value of the peak differ from the value of the 
plateau Xasikmin) only for a constant factor (see inset 
of Fig. 7) and therefore scales with the same exponent of 
the mean cluster size 7 [44] . Even if the long time limit 
of xj(a, t) is strongly different from the one observed in 
X4(fc,t), both fluctuations manifest a strong dependence 
on length scale. In fact, the peak of X?(fliO strongly 
decreases as a decreases (see Fig. 8). This may be inter- 
preted as a sign that clusters of bonded particle dominate 
the dynamics. 

Our data and these considerations suggest that het- 
erogeneities detected by x^(a,i) are due to the presence 
of clusters. However, despite the permanent nature of 
clusters, fluctuations of the self overlap x?('3-, t) decay to 
zero in the long time limit. This is due to the form of 
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FIG. 6: (Color online) Main frame: Self overlap, Qs{t), for 
a = 3 and different volume fractions from </!> = 0.02, 0.05, 
0.06, 0.07, 0.08, 0.09, 0.1 (from bottom to top). Inset: Self 
overlap for the same values of tj) of main frame and a = 0.15 . 



FIG. 7: (Color online) Main frame: Fluctuations 
of the self-overlap, x?('^i*)i for a = 3 and </> = 
0.02, 0.05, 0.06, 0.07, 0.08, 0.09, 0.1 (from bottom to top). In- 
set: 10-xJ(i*) (circles) for o = 3 and Xas(fcmin, 0) (triangles) 
as a function of {(be ~ </>). 



the overlap function u;(|fi(0) — ri{t)\), which is zero when 
a particle has moved a distance greater than a. There- 
fore two particles in the same cluster will contribute to 
X^(a,i), if the center of mass of the cluster has moved 
a distance less than a. In fact, when the cluster moves 
a distance larger than a, due to the form of the overlap 
function, the contribution to X?('^j ^) vanishes. Therefore 
we expect that for 1/a ~ kmm the peak is proportional to 
the mean cluster size, and occur at a time t* of the order 
of the time in which the center of the typical cluster of 
dimension ^ has moved a distance of the order of a. 

V. MEAN SQUARE DISPLACEMENT AND 
THE NON-GAUSSIAN PARAMETER 



Finally we have measured the mean square displace- 
ment (MSD) 

1 ^ 

^^'W = ]^E[(l^^'W-^^^(0)l')]' (12) 

i=l 

where fi{t) is the position of the i-th particle at the time 
t. In the main frame of Fig. 9 the MSD is shown for 
different volume fractions. Due to the Newtonian dy- 
namics, we find at very short time a ballistic behaviour, 
Ar^ (t) oc , followed by a crossover to a diffusive regime 
Ar^(t) oc t. The long time diffusive regime is always 
recovered for all the volume fractions considered, indi- 
cating that even at the percolation threshold this quan- 
tity is dominated by free motion of particles or clusters. 
Accordingly, no divergence of the inverse diffusion coef- 
ficient is found at the percolation threshold, where the 




FIG. 8: (Color online) Main frame: Fluctuations of the self- 
overlap, *)) for </> = 0.1 and a — 0.15, 0.5, 1, 1.5, 2, 2.5, 3 
(from bottom to top). 



numerous small size clusters continue to diffuse into the 
large mesh of the spanning cluster. 

We have also evaluated the MSD of the clusters and 
extracted their diffusion coefficient as a function of the 
cluster size s. In particular we obtained that for (f) — (pc 
D{s) for large s is fitted by a power law with 
/i = 1.0 ± 0.1 (see inset of Fig.9). Following [38] we 
expect D{s) - g^d-^+f /,.)/df ^ where d = 3 is the Eu- 
clidean dimension, / is the exponent which gives the di- 
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FIG. 9: (Color online) Main frame: Mean square displace- 
ment for cj} = 0.05, 0.06, 0.07, 0.08, 0.09, 0.1, 0.11, 0.12. 
Inset: Diffusion coefficient, D{s), as a function of the cluster 
size s for S = Sc. 



vergence of the viscosity, v ~ 0.88 is the critical exponent 
which gives the divergence of the connectedness length 
(see Appendix), and df ~ 2.4 is the fractal dimension 
of the spanning cluster at the threshold (see Appendix). 
Using these values we obtain a prediction for the expo- 
nent, which gives the divergence of the viscosity at the 
threshold / = v{hdf — d-f 2) ^ 1.23 in agreement within 
the errors with our data for the structural relaxation time 
(see Sect.III and Fig.2). 

In order to characterize the displacement of particles 
we have calculated the self part of the Van-Hove function 
[45]: 



G.(r,t) = l 



N 



{Y^5{r-\m-m\)) 



(13) 



If the motion of particles is diffusive with a diffusion co- 
efficient D, Gs{r,t) = (l/47ri:>t)3/2e(^'-V4i?t)^ where r is 
the distance traveled by a particle in a time t. In the in- 
set of Fig. 10 we plot the self van Hove function for a fixed 
volume fraction at different times. Our results indicate 
that for short times and short distances the function is 
well fitted by a Gaussian. For long distances and long 
times, the van Hove function is well fitted by an expo- 
nential decay. An exponential decay has been observed 
in different glassy systems for intermediate times [46]. 

The deviation from the Gaussian distribution at long 
times, observed in our system, indicates that some par- 
ticles move faster than others, due to the presence of 
clusters. Particles belonging to different clusters have 
a different diffusion coefhcient depending on the cluster 
size. As a consequence we suggest that, in the diffusive 
regime, Gs(r, t) does not have a Gaussian form, but it is 



FIG. 10: (Color online) Main frame: The self part of the 
Van-Hove distribution for = 0.09 and time t — 1285.02. 
The continuous line is obtained from the diffusion coefficient 
of clusters using Eq.(14). Inset: The self part of the Van- 
Hove distribution for ij) = 0.09 and different times t = 0.469, 
6.739, 93.199 (from left to right). The lines are Gaussian 
fitting functions. 



instead given by a superposition of Gaussians 



/ 1 



\ATTD{s)t 



3/2 



^-r-/4Dis)t^ (14) 



where D{s) is the diffusion coefficient of cluster of size s 
and ris is the cluster size distribution. 

In Fig. 10 we compare our data with Gs{r, t) calculated 
using Eq.(14) and D{s) obtained from the simulations. 
As we can see in figure, data are in good agreement with 
our predictions, provided that time is sufficiently long for 
clusters diffusing with diffusion coefficient D{s). 

In agreement with this finding, the non-Gaussian pa- 
rameter, which is a measure of the departure from the 
Gaussian behaviour of the probability distribution of the 
particle displacements, does not go to zero at long times. 
The non-Gaussian parameter is defined as [40] : 



a2{t) = 



5(Ar2(t))^ 



-1, 



(15) 



where Ar^{t) = ^ J2i^=i [Unit) ~ r\(0)|'')] and it is zero 
if the probability distribution of the particle displace- 
ments is Gaussian. 

In glassy systems [2] , (i) on the time scale at which the 
motion of the particles is ballistic, a2 is zero; ii) upon 
entering the time scale of the /? relaxation, a2 starts to 
increase; iii) on the time scale of the a relaxation, a2 de- 
creases to its long time limit, zero. The maximum value 
of a2 increases with decreasing temperature, signalling 
that the dynamics becomes more heterogeneous. 
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VI. CONCLUSION 
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FIG. 11: (Color online) Non-Gaussian parameter, Q2(t) as a 
function of time t ioi (j> = 0.05, 0.06, 0.07, 0.08, 0.09, 0.1, 
0.12. 



In the present model for permanent gels, we find that i) 
as in glasses, 012 is always zero on the time scale at which 
the motion of the particles is ballistic; ii) it tends in the 
long time limit to a plateau value, which increases with 
increasing volume fraction; iii) at low volume fraction, 02 
has a maximum at intermediate times, which disappears 
upon approaching the gelation threshold; iv) no critical 
behaviour is observed at the gelation threshold. 

Within our interpretation the asymptotic value of the 
non-Gaussian parameter, using Eq.(14) may be written 
in the following form: 



- 1 = 



D^-D 



(16) 



where, for each bond configuration, (...) is the average 
over the cluster distribution. We have verified that a2'' 

coincides with ^^—P within the errors. Hence our re- 

D 

suits indicate that the non-Gaussian parameter tends to 
a plateau given by the ratio of two quantities, which both 
have no critical behaviour at the percolation threshold. 
In summary, as the fluctuations of the self ISF, the non- 
Gaussian parameter does not decay to zero in the long 
time limit, due to the presence of permanent clusters. 
However, the main contribution to a2 comes from the 
numerous finite clusters (the bigger the cluster, the lower 
its diffusion coefficient D{s), and consequentially its con- 
tribution to 02), so that no critical! ty approaching the 
percolation threshold is observed in the non Gaussian 
parameter. 



We have presented a molecular dynamics study of 
a model for permanent gels and investigated its static 
and dynamical properties. Usually the sol-gel transition, 
marked by the divergence of viscosity and the onset of 
an elastic modulus, is interpreted in terms of the appere- 
ance of a percolating cluster of monomers linked by bonds 
[32, 33, 42]. While the viscosity and the elastic modu- 
lus can be measured directly, usually the experimental 
determination of percolative properties needs the manip- 
ulation of the sample (for a review see [47] and references 
therein), i.e. the sample is dissolved in a known quan- 
tity of solvent in such a way that each cluster is separated 
from the others. For the first time our results identify the 
thermodynamical observable associated with the cluster 
properties in a gelling system, and, via the measure of the 
fluctuations of the self ISF, allow to obtain the critical 
exponents without such a manipulation of the sample. 

In our model the formation of permanent bonds be- 
tween the particles leads to a percolation transition in 
the universality class of random percolation. The per- 
colation threshold coincides with the gelation threshold, 
marked by the slowing down of dynamics on length scale 
of the whole system. We have found that the behavior 
of the self ISF in the sol phase and near the threshold 
is in agreement with typical experiments on gelling sys- 
tems. In chemical gels the onset of a stretched exponen- 
tial decay is typically associated to the wide cluster size 
distribution close to the gelation threshold, producing a 
wide distribution of relaxation times. At the percolation 
threshold, the longest relaxation time diverges due to the 
critical growing of the percolation correlation length, pro- 
ducing a long time power law decay. Our results confirm 
this picture but new insights are obtained with a study 
of the dynamical heterogeneities, in terms of fluctuations 
of different correlation functions. In the present model 
for permanent gels, the fluctuations of the self-overlap 
present always a peak, whereas the fluctuations of the 
self ISF are monotonically increasing with time. Differ- 
ently from glassy systems, the fluctuations of the self ISF 
tend in the long time limit to a plateau, whose value, for 
the lowest wave vector, coincides with the mean cluster 
size. The behavior of the non-Gaussian parameter as 
a function of time is qualitatively similar: in the long 
time limit it reaches a plateau, due to the contribution 
of particles belonging to different clusters with a size de- 
pendent diffusion coefficient. Nevertheless, the value of 
the plateau does not diverge at the gelation transition, 
being dominated by the presence of small clusters with 
finite diffusivity. 

This study has shed some light on the differences be- 
tween the dynamics and the dynamical heterogeneities 
in glasses and chemical gels. We have been able to clar- 
ify that, when clusters of bonded particles are present, 
different time correlators can deliver very different infor- 
mation whereas in the studies on glasses they are often 
used interchangeably. On this basis, these findings have 
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FIG. 12: (Color online) Main frame: Percolation probability 
n ((/!>, L) as a function of the volume fraction tj) for boxes of size 
L = 15, 30, 40, 60. Inset: Data collapse obtained plotting 
U{<p,L) versus (0 - ^c)!/^^" with u = 0.88 and (pc = 0.09. 



FIG. 14: (Color online) Main frame: Average number of 
clusters (per particle) with masss versus sat for L — 40. 
The full line is a power law, s~^, withfitting parameter r = 
2.1. Inset: Radius of gyration, Rg, as a function of the mass 
s of clusters at <jic for L — 40. The full line is a power law, 
s^^'''f , with fitting parameter df = 2.4. 



^1000 
750 



^ 0.3 










^0.25 










0.2 

0.15 


' 'I 






L=15 


0.1 


V 






□ L=30 


0.05 








C L=40 




-0,2 0,2 


t 












1 


T L=60 




z □ 

, 1 , 5 , ^ , ■; , 


□ 


t 


h 

■,T , 



■•06 0.07 0.08 0.09 0.1 0.11 



FIG. 13: (Color online) Main frame: Mean cluster size 
S{(j}, L) as a function of the volume fraction (j) for boxes of the 
same sizes of Fig. 12. Inset: Data collapse obtained plotting 
S{(fi,L)L-^/'' versus (0 - (f)c)L^^'' with ly = 0.88, (f)c = 0.09 
and 7 = 1.85. 



der parameter" , whose fluctuations reveal the presence of 
heterogeneities in the dynamics, may be quite different. 
However both quantities show a strong length scale de- 
pendence (both strongly decrease as k or 1/a increase), 
which seems to be the distinct sign of (permanent or col- 
loidal) gelation compared with the (attractive or repul- 
sive) glass transition [31]. This result is confirmed by a 
recent work [48], where it has been found that, in a model 
for colloidal gels, at low volume fraction, the fluctuations 
of the self ISF for small wave vector display a depen- 
dence on time, which is dramatically different from the 
one found at higher volume fraction [23, 26] in the glassy 
regime. As a final remark, it is interesting to note that 
in the model here discussed the dynamical susceptibil- 
ity is similar to that observed in a spin glass model with 
quenched interactions [49] , suggesting a possible common 
description of the phase transition involved, as also pro- 
posed elsewhere [50]. 
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PCI and S.Co.P.E. 



interesting implications for the study of gels due to non- 
permanent bonds, as in the case of colloidal gels. In fact, 
our study also indicate a possible way to discriminate 
between a gel-like behaviour and a glass-like behaviour 
in these systems. Our results strongly suggest that, if 
heterogeneities are due to clusters of particles connected 
by permanent (or persistent) bonds, as in permanent (or 
colloidal) gels, the behavior of the "time-dependent or- 



APPENDIX A: PERCOLATION TRANSITION 

In this Appendix, with a finite size scaling analysis, 
the percolation threshold and the critical exponents are 
obtained. We find that the percolation of permanent 
bonds, corresponding to the sol-gel transition [32, 33], is 
in the universality class of random percolation. 

Varying the volume fraction 0, we have measured the 
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percolation probability, 11(0) (defined as the average 
number of configurations where a percolating cluster is 
found), the cluster size distribution, n,,, and the mean 
cluster size S{(j)) = ^s^ns/^sng. For each volume 
fraction we have used simulation boxes of different size L 
and, from a standard finite size scaling analysis [42], we 
have obtained the percolation threshold (/>c, and the criti- 
cal exponents v (which governs the power law divergence 
of the connectedness length ^ ^ — (pc]^'^ the tran- 
sition threshold is approached from below) and 7 (which 
governs the power law divergence of the mean cluster size 
S ^ \(j) — (j)c\~'')- The percolation threshold and the criti- 
cal exponents obtained from the data showed in Fig.s 12 
and 13, are respectively <j)c 0.09 ±0.01, p = 0.88 ±0.05 
and 7 = 1.85 ± 0.05. The cluster size distribution 



for (f) = (j)c, shown in main frame of Fig. 14, follows a 
power law behavior ~ with a Fisher exponent 
r = 2.1 ±0.2. 

Finally in the inset of Fig. 14 the radius of gyration 
Rg as a function of the mass s of clusters is plotted. 
The data are well fitted by a power law with expo- 
nent 1/df = 0.42 ± 0.03 (full line in figure) which gives 
df = 2.4 ± 0.1 in agreement with the fractal dimen- 
sion of the random percolation clusters in Sd, df ~ 2.5. 
The measured values of the critical exponents satisfy the 
hyper-scaling relations (2/3 + 7 = i^d, dj = d — P/v, and 
T = 2 + {d — df)/df [42]), and are in good agreement with 
those of the 3d random percolation {v = 0.88, 7 = 1.80 
and T = 2.18 [42]). 
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